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ABSTRACT 

Infrared extinction maps and submillimeter dust continuum maps are pow- 
erful probes of the density structure in the envelope of star-forming cores. We 
make a direct comparison between infrared and submillimeter dust continuum 
observations of the low-mass Class core, B335, to constrain the ratio of submil- 
limeter to infrared opacity (n S mm/ Kir) an d the submillimeter opacity power-law 
index (k oc A - ^). Using the average value of theoretical dust opacity mod- 
els at 2.2 /iin, we constrain the dust opacity at 850 and 450 fim. Using new 
dust continuum models based upon the broken power-law density structure de- 
rived from interferometric observations of B335 and the infall model derived 
from molecular line observations of B335, we find that the opacity ratios are 
^ = (3.21 - 4.80) l^; 44 , x 10" 4 and ^ = (12.8 - 24.8)± 2 ; 4 x 10~ 4 with a submil- 
limeter opacity power-law index of (3 smm = (2.18— 2.58)1q;3q. The range of quoted 
values are determined from the uncertainty in the physical model for B335. For 
an average 2.2 /zm opacity of 3800 ± 700 cm 2 g _1 , we find a dust opacity at 850 
and 450 /mi of k 850 = (1.18 - 1.77)1^ and k 45 o = (4.72 - 9.13)±J| 8 cmV 1 of 
dust. These opacities are from (0.65 — 0.97)k° 5 q 5 of the widely used theoretical 
opacities of Ossenkopf and Henning for coagulated ice grains with thin mantles 
at 850/im. 
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1. Introduction 

The commissioning over a decade ago of two dimensional bolometer cameras such as 
SCUBA (Submillimeter Common User Bolometer Array; Holland et al. 1999) and SHARC 
(Submillimeter High Angular Resolution Camera; Hunter et al. 1996) permitted the efficient 
mapping of dust continuum emission at submillimeter wavelengths. Submillimeter contin- 
uum observations of star forming regions in the Milky Way have constrained the physical 
structure of dense, star-forming cores through all embedded phases of core and protostar for- 
mation (e.g., Shirley et al. 2002; J0rgensen et al. 2002; Williams, Fuller, & Sridharan 2005, 
Kirk et al. 2005). On larger scales, mapping of entire molecular clouds has constrained the 
the dense core Initial Mass Function (IMF) and made surprising connections to the shape of 
the stellar IMF (e.g., Motte et al. 1998, Johnstone et al. 2001, Enoch et al. 2007, Motte et al. 
2007). In all of these studies, the mass or density scale is set by the assumed submillimeter 
dust opacity k (cm 2 g _1 ) since the mass of optically thin submillimeter emission is inversely 
proportional to the dust opacity (M d oc 1/k; Hildebrand 1983). It is very important to use 
an accurate value of the dust opacity since the mass distribution within the protostellar core 
directly affects the dynamical stability of the core as well as the radiative transfer through 
the core. The uncertainty in submillimeter dust opacity is the largest source of uncertainty 
in mass calculations and radiative transfer models of protostellar cores (see Shirley et al. 
2005). 

Since the pioneering work of Knacke & Thompson (1973), there have been many at- 
tempts to estimate the submillimeter dust opacity in star forming cores. The opacity of 
dust grains in the general ISM can vary substantially from dust opacities in environments 
surrounding star formation. From optical studies of dust absorption, it is well known that 
the ratio of total to selective extinction, Ry = Ay/E(B — V), varies from Ry = 3.1 in the 
general ISM to larger values (Ry = 5.5) in denser star forming regions (e.g., Mathis 1990, 
Whittet 2003, Draine 2003). Several important physical processes directly affect the opaci- 
ties in dense regions (Henning, Michel, & Stognienko 1995). The compositions of dust grains 
may vary from region to region. Dust grains may coagulate, changing the shape and the 
normalization of the general ISM size distribution (e.g., Mathis, Rumpl, & Nordsieck 1977; 
Mathis, Mezger, & Panagia 1983, Ossenkopf 1993, Ormel et al. 2009). Most of the molecular 
gas in a dense core is shielded from energetic photons from the forming protostar or from 
the Interstellar Radiation Field (ISRF) and confined to temperatures less than 20K. Many 
molecular species freeze out of the gas phase forming layers of polar (H 2 0) and apolar (CO) 
ices that change the dielectric properties of the grains and the size distribution of grains. 
Theoretical calculations of the opacities make various assumptions about grain composition, 
grain size distributions, grain geometry, grain porosity, and ice mantle compositions; the 
resulting predicted submillimeter opacities vary by up to an order of magnitude (see Table 
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2 of Shirley et al. 2005). Currently, the most widely used calculation of opacities is that of 
Ossenkopf & Henning (1994) which takes into account coagulation and varying thicknesses 
of ice mantles for dust grains that have persisted at high density, n = 10 6 cm -3 , for 10 5 years. 
Even the Ossenkopf & Henning opacities vary by a factor of a few at submillimeter wave- 
lengths depending on the particular assumptions used in the model and there is no guarantee 
that these opacities are appropriate for the wide variety of environments towards which they 
have been applied, from cold, low-mass starless cores to warm, high-mass proto-cluster cores. 

Clearly, observational constraints on submillimeter dust opacities are needed. There 
have been a few attempts to constrain the dust opacity at long wavelengths. These meth- 
ods use observations of the amount of dust extinction at near-infrared wavelengths directly 
compared to the amount of dust emission at (sub) millimeter wavelengths to constrain the 
opacity ratio between submillimeter and near-infrared (or visible) wavelengths. Kramer et 
al. (1998,2003) studied the opacity ratio in several dense starless cores in the IC 5146 region. 
Similarly, Bianchi et al. (2002) studied the nearby low-mass starless core B68. Both groups 
have demonstrated that the technique works if sensitive observations are obtained at both 
submillimeter and near-IR wavelengths and care is taken to compare observations taken at 
two very different resolutions (pencil beam at near-infrared wavelengths versus beam con- 
volved emission at submillimeter wavelengths with non-Gaussian telescope power patterns). 
This method has never been applied before to a dense core which harbors a protostar 

The near-infrared to submillimeter comparison requires high signal-to-noise observations 
at three wavelengths (two near-infrared wavelengths and one submillimeter wavelength). 
Few submillimeter studies of protostellar regions have sensitive maps that detect extended 
emission at high signal-to- noise ratios over large regions of the core (> 2') because the 
ground-based observations are limited by the size of the chopping needed for sky subtraction. 
One promising object is the nearby, well studied Class protostar located within the dense 
isolated Bok globule, Barnard 335. Since its initial detection as a far-infrared source (Keene 
et al. 1983), B335 has received considerable attention as it is one of the best protostellar 
infall candidates (Zhou et al. 1993, Choi et al. 1995, Evans et al. 2005) identified during 
the deeply embedded phase of low-mass star formation. For the purposes of this study, B335 
has been observed with high sensitivity with SCUBA at submillimeter wavelengths (Shirley 
et al. 2000) and with NICMOS at near-infrared wavelengths (Harvey et al. 2001). 

In this paper, we utilize a method of comparing the near-infrared extinction to submil- 
limeter emission to constrain the dust opacity ratio between submillimeter wavelengths and 
2.2 fim (§2) toward B335. From dust continuum radiative transfer we derive an updated 
physical model for B335 using currently published dust opacities (§3). The dust opacity 
ratio is determined utilizing our constraints on T^(r) (§4.1). We then constrain the opacity 
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at submillimeter wavelengths and the opacity power law index, f3 smm (§4.2,4.3). Throughout 
this paper, we refer to the dust mass opacity, k u , in units of cm 2 per gram of dust. 



2. Method 

By directly comparing infrared extinction to submillimeter emission maps along com- 
mon lines of sight, we constrain the opacity ratio K srnrn /K ir . The emission at submillimeter 
wavelengths along a line-of-sight, s, is given by the equation of radiative transfer for optically 
thin emission, dl v jdr v = B v [T(s)] with 

s ^ m =nm H (m d /m g )K smm (s) n(s) , (1) 

where B u is the Planck function and n(s) is the gas number density (cm~ 3 ). We use the 
dust mass opacity, k v (cm 2 g _1 ), where we have assumed a mean molecular weight \x = 2.32 
and an average gas mass to dust mass ratio of (rrid/m g ) = 1 : 100. In the derivation of 
the expression to determine the opacity ratio, \i and {md/m g ) cancel out; but they are 
stated explicitly here since they are used in radiative transfer modeling of the submillimeter 
emission (§3). Integrating along the line-of-sight gives 

I smm 

(s)ds , (2) 

(Adams 1991, Shirley et al. 2000). 

At infrared wavelengths, the observed intensity is due to the total amount of extinction 
along the line-of-sight, 

I ir (s) = 4-(0)e-^ ; (3) 
therefore, the total extinction in magnitudes is 

A ir = -2.5 log (yY^j^j = 2.51og(e) J^m H (m d /m g )n(s)K ir (s)ds . (4) 

The infrared opacity includes contributions from absorption and scattering (/Cj r = k,^ s + K,^ a ) ■ 
Dividing Equation (2) by Equation (4), we derive the relationship between the submillimeter 
intensity, the infrared extinction, and the opacity ratio 

T _ j s B u [T{s)}n{s)K Sill III {s)ds 

lsmm ~ 2.5\og(e) J s n(s)K tr (s)ds Ar " ^ 

If we further assume that the infrared and submillimeter opacity does not vary along the 
line-of-sight, then we find 

lsmm ( K>smm \ , . 

A ir , (6) 
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where P n is related to the density-weighted average Planck function and is defined as 

f s n(s)B u [T(s)]ds 
2.5 log(e) J n(s) ds 

If an isothermal approximation is used, then the column density cancels in P n and we are 
left with P n = f?„(Tj SO )/2.51og(e). However, dust temperature gradients exist throughout 
protostellar envelopes (Shirley et al. 2003). In this paper, we shall use the n(r) and T<j(r) 
determined from the best-fitted ID dust continuum radiative transfer models to calculate 
P n along each stellar line-of-sight. 

Theoretically, the opacity ratio is determined from the slope of a plot of the submil- 
limeter intensity versus the near-infrared extinction. In reality, neither the submillimeter 
intensity nor the near-infrared extinction are directly observed. Submillimeter observations 
actually observe the convolution of the source specific intensity distribution with the tele- 
scope beam pattern such that the observed flux density is S smm = (I smm )£lbeam- A typical 
submillimeter beam pattern is not well described with a single Gaussian main beam as a 
significant fraction of the power pattern is contained within the sidelobes. The total solid 
angle of the beam, f2 beam , is determined from the integral of the normalized telescope power 
pattern, including sidelobes, over solid angle. 

With observations of two or more near-infrared wavelengths, we observe color differences 
of background stars. The extinction along a line of sight is determined from a scaling law 
that relates A{ r to the observed color excess. In this paper, we will determine the extinction 
at the near-infrared K band at 2.2 /im. We define the K band selective extinction, R K such 
that 

A 2 . 2 = Rk[(H -K)- ((H - K) }} = R K E(H - K) = R K (A 1S5 - A 2 . 2 ) , (8) 

where ((if — K) ) is the mean intrinsic (H-K) color of the background stars and (H — K) 
is the observed infrared colors with extinction. The determination of Rk assumes that the 
near-infrared extinction law is well described by a power-law (A{ r oc \~P™ r for near-infrared 
wavelengths; see Draine 2003 and references therein, also Flaherty et al. 2007). Chapman & 
Mundy (2009) determine (3 nir = 1.7 from the slope of a plot of E(J — H) versus E(H — K) of 
background stars toward four dense cores. This is in the middle of the range of the typical 
values of /3 n j r = 1.6 — 1.8 found in the literature (Draine 1989, Rieke & Lebofsky 1985, 
Martin & Whittet 1990, Whittet et al. 1993). Using this range of (5 nir in Equation (8), we 
find R K = 1.59 ± 0.12. We shall use this value of R K throughout this paper. 

We constrain K smm / ' k 2 . 2 from the slope of a plot of S smm / P n VL heam versus E(H — K). 
If the slope is denoted by b, then the opacity ratio is simply n smm /K 2 , 2 = b/Rx- It is 
important to reiterate that this derivation assumes that the dust opacity ratio along an 
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individual line-of-sight is constant. If there is a variation in the dust opacity, then a linear 
regression is measuring an emission-weighted average dust opacity ratio along each line-of- 
sight. A monotonic change in the opacity ratio with radius will produce curvature in a plot of 
Ssmml Pn^beam versus E(H — K) while a distribution of opacity ratios along each line-of-sight 
will produce intrinsic scatter in the correlation. 

2.1. Submillimeter Images 

The reduction and analysis of SCUBA 850 and 450 /im jiggle maps is described in detail 
by Shirley et al. (2000). We have re-analyzed and combined SCUBA images of B335 taken 
by Shirley et al. (2000) and images from the SCUBA CADC archive. B335 was observed 
with SCUBA in jiggle mapping mode on only two nights (April 17, 1997 and December 18, 
1997) with low atmospheric opacity (r 2 25 < 0.05). Throughout this paper, for simplicity we 
quote the SCUBA wavelengths as 850 and 450 /im, whereas the actual narrow-band SCUBA 
filters have average wavelengths of 860 and 445 /zm respectively. The average wavelengths 
are not sensitive to the shape of the source spectrum (i.e. the average wavelength changes 
from 860 to 859 /im for sources with u° to z/ 4 .) 

Before combining images from two different nights, observations of Uranus taken within 
1 hour before and after the B335 observations on each night were analyzed to compare 
the shape of the telescope beam pattern. We found no significant difference in the Uranus 
radial profiles between the two sessions. The flux density scale is calibrated using the peak 
and integrated Uranus flux observed on the same night as the core as determined from the 
FLUXES program % We found very good agreement in the peak flux on both nights. The 
specific intensities in the final image (erg s _1 cm -2 Hz -1 sr -1 ) were estimated from the 
measured flux densities by using the average of the solid angle of telescope beam pattern 
determined on the same nights as the B335 observations ((f^sso) = 7.3 x 10 -9 sr and (fi 450 ) = 
4.4 x 10 -9 sr). The rms noise in the combined images is 15 mJy/beam at 850 /im and 82 
mJy/beam at 450 /im. The combined images are very similar to the published contour maps 
in Shirley et al. (2000; see Figure 1). 

The final image pixels are over-sampled (1") to determine the closest lines-of-sight for 
comparison with background stars (Figure 1). Since multiple background stars may lie within 
a single SCUBA beam (15" at 850 /im and 8" at 450 /im), the submillimeter intensity will be 
semi-correlated in a plot of S smm /P n fl versus E(H — K). By using the oversampled SCUBA 
map, we preserve all of the information in the near-infrared extinction map. Smoothing 
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the near-infrared line-of-sight to a 20" — 30" resolution, as is typically done with extinction 
mapping methods (e.g. Teixeira et al. 2005; Alves et al. 2001), results in a suppression 
of the true scatter in line-of-sight E(H — K). For this reason, we do not smooth the near- 
infrared data. Similarly, attempting to deconvolve the SCUBA beam from the B335 image 
is extremely difficult because the two dimensional beam shape varies during the observations 
due to changes in the shape of the telescope surface and variation of parallactic angle. 

2.2. Near-infrared Images 

The near-infrared images are from NICMOS observations and are analyzed in detail by 
Harvey et al. (2001). Observations were performed with the NIC3 camera on the Hubble 
Space Telescope using the F160W and F222M filters. A detailed photometric comparison 
was made between the NICMOS magnitudes and J, H, and K observations made with the 
NIRC camera on Keck I. B335 was imaged with a 3 x 3 mosaic (NIC3 FOV 51"2) plus a 4' 
radial strip centered on the Class protostar. Background fields, off the Barnard 335 cloud, 
at similar galactic latitude were observed with NICMOS and NIRC. The mean intrinsic 
color of background stars is ((H — K) ) = 0.13 mag and intrinsic scatter in the color is 
<j{H-K) = 0.16 mag. 

The distribution of background stars compared to 850 /iin emission is shown in Figure 
1. More than 200 background stars are observed. Most of the background stars are located 
outside of 20" with only 4 stars observed within 20" of the protostar and none within 15". 
The inability to see background stars in the innermost regions of Barnard 335 is due to the 
high column densities observed toward the core. Therefore, a submillimeter to near-infrared 
opacity comparison can only be made in the outer regions of the envelope greater than 15" 
from the protostar (3750 AU projected at a distance of 250 pc). 

3. An Updated ID Dust Model for B335 

The density and dust temperature along each line-of-sight are determined from radiative 
transfer models of the dust continuum emission (e.g., Shirley et al. 2003). The radiative 
transfer models self-consistently calculate the dust temperature profile, T^(r), using a one 
dimensional radiative transfer code (CSDUST3; Egan, Leung, & Spagna 1988) determined 
from an input density distribution (n(r)), interstellar radiation field, internal luminosity, and 
a dust opacity curve The model intensity profiles and spectral energy distribution 

are reconstructed using the same techniques as the observations (e.g., beam convolution, 
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chopping, aperture matching). The best-fitted models minimize the f° r the observed 
submillimeter intensity profiles and the observed spectral energy distribution at wavelengths 
where the optical depth is less than unity (A > 60 /mi). Details of the radiative transfer 
modeling procedures for low-mass cores may be found in Shirley et al. (2002, 2005). 

Unfortunately, in order to calculate P n and ultimately the observed opacity ratio, we 
have to assume a dust opacity curve to input into the radiative transfer model. Therefore, 
we must explore different theoretical opacities to determine how our derived opacity ratio 
is biased by our radiative transfer opacity choice. Since the background stars toward B335 
are at projected lines-of-sight greater than 1000 AU from the protostar and since B335 is 
a low luminosity protostar (L^ = 3.3 L ) embedded in a dense core that is exposed to a 
weaker than average ISRF (Shirley et al. 2002; Evans et al. 2005), the calculated T(r) at 
those radii will be a slowly varying quantity with radius (Shirley et al. 2002). In this paper, 
we explore the effect on T^(r) for different dust opacity models calculated by Ossenkopf & 
Henning for coagulated grains with varying thicknesses of ice mantles (columns 2,5, and 8 of 
Ossenkopf & Henning 1994) and dust opacity models calculated by Weingartner & Draine 
(2001) for ISM grain populations with R v = 3.1, 4.0, and 5.5 (size distribution "A"@). 
Since the original Ossenkopf & Henning models did not calculate the scattering opacity, 
the OH opacity was divided between scattering and absorption using the ratios from the 
Pollack et al. (1994) models that best match the Ossenkopf & Henning absorption opacity 
(Young & Evans 2005). The ratio of scattering to absorption opacities across the 3 /jm ice 
feature were determined from the albedos in Pendelton et al. (1990, Figure 4b) since the 
Pollack et al. models did not include ices. A complete explanation of the modifications to 
the short wavelength OH opacities may be found in §2.1 of Young & Evans (2005). These 
short wavelength-modified OH opacities have been used in several published dust continuum 
radiative transfer calculations (e.g., Young & Evans 2005, Shirley et al. 2005, Dunham et 
al. 2006, Dunham et al. 2010). 

We assume a distance of 250 pc (Tomita et al. 1979) to be consistent with previously 
published models, although this distance is very uncertain (see Olofsson & Olofsson 2009). 
Shirley et al. (2002) explored the effects of a closer distance (125 pc) on the radiative transfer 
models; however, the determination of the opacity ratio does not depend on distance. 

There is also uncertainty in the modeled density structure of the B335 core. The original 
published dust continuum models of B335 found that a single power-law density model 
(n(r) oc r~ L8 ) using OH5 dust opacities is a good fit to the observed intensity profile while 
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the fit to the shape of the far- infrared SED was not as well matched (Shirley et al. 2002). 
Subsequently, Harvey et al. (2003a,b) published an updated model based on interferometric 
millimeter continuum emission and near-infrared extinction maps (Harvey et al. 2001) and 
found that a broken power-law was a better fit to the density structure with a flatter density 
profile (r~ L5 ) inside 6500 AU and steeper profile outside (r~ 2 '°). Most recently, Doty et 
al. (2010) used broken power-laws and variable opacities in four radial zones to attempt to 
model the large-scale dust continuum emission (excluding interferometric constraints). All 
of these modeling efforts found strong disagreement with the previously published molecular 
line radiative transfer modeling which is best-fit by a Shu (1977) infall solution with a 
modest infall radius r in f = 6200 AU (Choi et al. 1995). Evans et al. (2005) confirmed this 
disagreement with updated radiative transfer modeling using non-uniform abundance profiles 
for several molecular species. Reconciling the modeling differences requires simultaneous 
modeling of molecular line and dust emission (interferometric and single-dish) with varying 
dust opacities with radius (e.g., including CO desorption for Td > 20 K). Unfortunately, this 
is beyond the scope of this paper; therefore, we shall analyze the uncertainty on the opacity 
ratio due to our uncertainty in the underlying density model by using both the best-fit 
molecular line model and the best-fit dust continuum model. 

We explore a grid of models that varies the dust opacity and the density scale factor, /. 
The entire density profile is scaled by a single number, /, to match the observed flux in a 120" 
aperture at 850 fim. The strength of the ISRF is constrained from molecular line modeling 
of CO observations (Evans et al. 2005). We adopt an ISRF parametrized in Shirley et al. 
(2005) and corresponding to G = 0.1 Habings (s isr f = 0.3, A V (R ) = 1.0 mag) at an outer 
radius of R D = 3 x 10 4 AU, consistent with the extent of the near-infrared extinction profile 
(Harvey et al. 2001). Two physical models are used: a scaled Shu (1977) infall solution 
(n(r) = /nshu( r )) with r in f = 6200 AU; and a broken power-law solution derived by Harvey 
et al. (2003a,b), 



Results for a subset of the models are summarized in Table 1. 

As noted by Shirley et al. (2002), the scaled Shu infall models do not fit the observed 
submillimeter intensity profiles or the SED. In contrast, the broken power-law of Harvey et 
al. scaled in density by a factor of / = 2.4 with OH8 opacities provide a good fit to the 
submillimeter intensity profiles and bolometric luminosity and a slightly better fit to the 
shape of the SED than the originally published best-fitted model by Shirley et al. (2002; 
Figure 2) which use OH5 opacities. All of the radiative transfer models have difficulty fitting 




r G [100,6500] AU 



r e [6500, 30000] AU . 



(10) 
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the far-infrared SED indicating that there may be a problem with the theoretical opacities 
at those wavelengths or effects of non-spherical geometry for B335. Since the Weingartner 
& Draine opacities are much smaller at submillimeter wavelengths than the Ossenkopf & 
Henning opacities, the density scaling factor is higher for WD models. We calculate P n (0) 
for all of the models listed in Table 1 to analyze the effect of the model opacities on the 
derived opacity ratio (§4.2). 

4. Results 

4.1. Linear Regression Technique 

In order to determine the opacity ratio between a submillimeter wavelength and 2.2 /zm, 
we must determine the slope from the plot of submillimeter intensity versus the near-infrared 
color excess. We use the Bayesian linear regression routine LINMIX_ERR (Kelly 2007) to 
determine the slope of a relationship of the form 

y = a + bx + a 2 int (11) 

with intrinsic scatter about the line, af nt , and heteroscedastic errors in both x and y. In our 
analysis, x = E(H — K) and y = S smm /P n Qbeam- LINMIX_ERR approximates the distri- 
bution of the independent variable ( mixture of Gaussians. This method alleviates 
the ad hoc assumption of a uniform prior distribution on the independent variable that is 
used in the derivation of popular x 2 minimization routines such as XYEFIT (e.g., Press et 
al. 1992, Tremaine et al. 2002, Weiner et al. 2006) and also permits fitting of truncated 
data sets (e.g., Malmquist bias) and data sets that include censored data or upper limits. 
Details of the assumed prior distributions are described in detail in Kelly (2007). Direct 
computation of the posterior distribution is too computationally intensive; therefore, ran- 
dom draws from the posterior distribution are obtained using a Markov Chain Monte Carlo 
method (Metropolis-Hastings Algorithm, Chib & Greenberg 1995; Tierney 1998; Gelman et 
al. 2004). We fit our data using the publicly available IDL code LINMIX_ERR.pro (Kelly 
2007) to determine the distributions of a, b, and a\ nt . In order to test the robustness of 
the resulting distributions, we varied the number of Gaussians from k=2 to k=4 and used 
various numbers of iterations up to 10 4 . 

The plots of S smm / P n Qb eam versus E(H — K) are shown in Figures 3 and 4 for the best- 
fitted dust continuum model (2.4 x Harvey BPL OH8) at 850/im and 450/im. The histograms 
of the slope and intercept are well approximated by Gaussian distributions (see inset of 
Figures 3 and 4); therefore, we tabulate the mean and standard deviations of the slope and 
intercept distributions (Table 1). The linear regression is only performed on data points 
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that lie outside the outflow cavities. Harvey et al. (2001) noticed a bimodal distribution 
of color excesses depending on whether a background star was located in the direction of 
the east-west oriented molecular outflow cavity. The opening angle of the CO outflow was 
recently characterized by Stutz et al. (2008) to be 55°. This is slightly larger than the 
outflow opening angle of 40 ± 5° used in Harvey et al. (2001) and (2003a). Since we are 
making a comparison of the opacity properties of the dust in the envelope, and since shock 
processing of dust within the outflow cavity walls may affect the grain opacities, we use the 
larger estimate for the outflow opening angle (55°) in this paper. 



4.2. Opacity Ratios and Submillimeter Opacity 

The opacity ratio is determined from plots similar to Figures 3 and 4 for all of the models 
discussed in §3. Linear regression from the LINMIX_ERR method provide good fits to the 
observed correlations. We find no evidence for large, systematic curvature in the plots of 
S smm /P n Q versus E(H — K) indicating that over the range of impact parameters probed by 
background stars (9 e [15", 70"]), there is not a strong monotonic gradient in the intensity- 
weighted dust opacity with radius. This result independently confirms the conclusions by 
Shirley et al. (2002) that there is no evidence for large scale opacity changes from comparisons 
of 850 and 450 /im model intensity profiles. This result is also consistent with the findings of 
Doty et al. (2010) which indicate nearly constant dust opacities over the range of radii that 
are fit in this paper. A caveat is that there is significant scatter in the correlations at both 
submillimeter wavelengths that could mask opacity changes. The scatter appears to become 
slightly larger near the highest E(H — K) observed. This may indicate that opacities are 
beginning to change within a few thousand AU of the protostar. Theoretically, we expect 
changes in the ice mantle composition as desorption due to protostellar heating occurs in the 
inner envelope. How this affects the observed dust opacity is still unknown. Unfortunately, 
NICMOS was not sensitive enough to detect to detect background stars within 10" of the 
protostar. 

The underlying density and temperature distribution affect the slope of the correlation 
through P n . In order to analyze the effect of the temperature profile on the opacity ratio, 
we first assume an isothermal approximation where every line-of-sight is assumed to have 
the same dust temperature (T(r) = constant). The calculated opacity ratios for isothermal 
temperatures are shown as the solid curve in Figure 5. For example, the 850 /im ratio 
varies by a factor of 2.5 for dust temperatures from 7 to 14 K. In reality, each line-of-sight 
through the protostellar envelope is non-isothermal (§3). We calculate P n (9) from Equation 
7 using the density (n(r)) and dust temperature profiles (T^(r)) determined from the dust 
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radiative transfer models for each line-of-sight (6). For each of the dust continuum models, 
we estimated the isothermal temperature through the stellar lines-of-sight used in the linear 
regression by solving Equation 7 for the temperature and using the calculated line-of-sight 
P n {9) from the radiative transfer model 



This temperature corresponds to the single temperature that characterizes a non-isothermal 
line-of-sight. The average (TJ os ) is calculated by averaging T/ os for the 190 lines-of-sight used 
in the linear regression (§4.1). When the more realistic T{r) from the radiative transfer 
models is included in the regression, then the calculated opacity ratio is always below the 
opacity ratio determined from isothermal lines-of-sight (Figure 5). This is a systematic 
effect caused by a monotonically decreasing temperature profile (dT/dr < 0; see Figure 2a) 
along each line-of-sight. It is very important to account for temperature gradients when 
determining the opacity ratio in the envelopes of Class protostars. 

Figure 5 also graphically illustrates the range of uncertainty introduced into the determi- 
nation of K smm / K ir due to different model opacity assumptions and different physical models. 
In quoting our opacity ratio, we choose two limiting models that characterize the range of 
K S mm/ and which also fit the observed 850 /im flux. Those two models are 2.8 x Harvey 
BPL OH5 and 2.8xEvans Shu OH8 in Table 1. We find that = (3.21-4.80)1^ x 10~ 4 
and = (12.8 — 24.8) l^g x 10~ 4 , where the range in K smm / n ir corresponds to the value 
for each physical model. 

We may compare our value of the opacity ratio with previous determinations toward 
dense cores (Figure 6). Bianchi et al. (2003) observed the starless core B68 at 850 and 1200 
/im with SCUBA and SIMBA. They employ a similar technique to compare the submil- 
limeter intensity and near-infrared colors. However, they assume isothermality of the dust 
temperature and the extinction law of Rieke & Lebofsky (1986). Since they determine the 
opacity ratio with respect to the opacity at optical wavelength, V, we must use the Rieke 
& Lebofsky A K /A V = 1/8.9 to convert from k v to n 2 .2- The resulting opacity ratios are 
^850/^2.2 = 3.6±0.9x 10" 4 and «i20o/«2.2 = 8.0±2.7x 10~ 5 . The 850 /zm ratio is comparable 
to the ratio we determined for B335. 

Similarly, Kramer et al. (2003) determined the 850 /xm opacity ratio toward 4 cores in 
the IC 5146 filament. This study includes an analysis of dust temperature variations between 
the cores in the filament. Again we must use the Rieke & Lebofsky A K /A V to convert to 
K2.2- The four cores have k^q/k 2 .2 that range from 1.9 ± 0.2 x 10~ 4 to 5.4 ± 0.3 x 10~ 4 . 
Kramer et al. also find evidence that the opacity ratio has an inverse dependence on the 




hu/k 
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dust temperature (see Figure 5). The (^850/^2.2) from all three studies is 3.7 ± 1.4 x 10~ 4 . 
Unfortunately, no opacity ratios have been determined previously in the literature at 450 
/im toward low-mass dense cores. 

Converting the observed opacity ratio to the submillimeter opacity requires an estimate 
of the opacity at 2.2 /im. A first approximation is to average the H2.2 from many different 
theoretical opacity models. Averaging the 2.2 /im opacities from the Ossenkopf & Henning 
models (OH2, OH5, OH8; 1994), the opacities used in the multi-dimensional dust models 
of Whitney et al. (private communication; e.g., Whitney et al. 2003), the Mathis et al. 
opacity (1983), and the new theoretical opacities calculated by K. Pontoppidan (private 
communication) that match the Cores to Disk mid-infrared extinction law and ice features 
(Pontoppidan et al., in preparation), we find (^2.2) = 3800±700 cm 2 g -1 of dust. Multiplying 
this number into the opacity ratios results in the submillimeter opacities of k 850 = (1.18 — 
1.77) ±0:24 and K450 = (4.72 - 9.13)lJ;g 8 cm 2 g" 1 of dust. These opacities are plotted with 
theoretical curves in Figure 6. We note that this crude average for K2.2 results in a large 
errorbar in the calculated 850 and 450 /im opacities because the 2.2 /im opacities vary by a 
factor of two among the different theoretical models. 

The uncertainty in the 850 and 450 /im opacity ratios and opacities make our de- 
terminations consistent with the empirical opacity law parameterized by Mathis (1990; 
k = 13.16(A/250/mi)~ 2 ; also parameterized in Kramer et al. 2003) as well as the coagu- 
lated dust model with thin ice mantles of Ossenkopf and Henning (OH5, 1994). The opacity 
model, OH8, that provides the best fit to the submillimeter and intensity profile and SED 
is at the upper statistical errorbar at 850 /im and the lower bound of the models at 450 /im. 
Our results bracket the popular theoretical opacities (OH5) at 850 /im that have been used 
in dust continuum radiative transfer modeling (e.g., Shirley et al. 2002, Mueller et al. 2002, 
Young et al. 2003, Shirley et al. 2005, Dunham et al. 2006, Doty et al. 2010). 

4.3. Power-law index f3 

At far-infrared and submillimeter wavelengths > 100 /im, the opacity falls as a power- 
law with increasing wavelength (k(A) oc A _|8smm ). Estimating /3 smm is a difficult problem. 
The most traditional methods have used modified blackbody fits to the SED or ratios of sub- 
millmeter wavelengths to constrain /3 smm (e.g. Visser et al. 1998, Shirley et al. 2000). Both 
of these methods assume a single dust temperature which is not an appropriate assumption 
for Class protostars which have strong temperature gradients (i.e., Fig 2). Instead, we 
use the derived opacity ratios at two submillimeter wavelengths to constrain f3 smm between 
450 and 850 /im. Since this method utilizes background stars that are at least 15" from 
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the central protostar, the opacity ratios are probing the dust properties in the outer, cold 
portions of the envelope. This differs from the previous two methods which use fluxes that 
include significant contributions from warmer dust near the protostar. The opacity ratio 
/3 smm is given by 



Since this is a ratio, the exact value of Rk cancels; however, the ratio is still sensitive 
to systematic uncertainties such as flux calibration errors. The opacity ratio for B335 is 



If we compare our results to Bianchi et al. (2003) which determined the opacity ratio 
at 850 and 1200 /im, we find a severe discrepancy that illustrates the importance of the 
calibration in determining (3. The Bianchi [3 srnrn = 4.3 ± 1.3 is much higher than our (3 smm . 
This anomalous result indicates a systematic calibration problem at one or both wavelengths. 
Great care must be taken when comparing observations made with different instruments 
on different telescopes through different observing conditions (e.g. SCUBA and SIMBA). 
Observations taken with SCUBA simultaneously at 850 and 450 /xm avoid this problem since 
the observations are taken in the same atmospheric conditions. While a calibration error 
at 450 or 850 /xm could account for our f3 smm > 2, we have taken great care to assure a 
stable calibration between 850 and 450 /im by comparing the flux calibration of Uranus 
taken on several nights surrounding the B335 observations. The ratio of the 450 to 850 /im 
flux conversion factors (see Jenness et al. 2002) never vary by more than 10% during these 
time periods. 

Observational evidence for submillimeter opacity indices above two in the ISM exist. 
A detailed multi- wavelength study of the starless core, TMC-1C, using multiple methods 
to determine the opacity index find that 1.7 < (3 < 2.7 with a most likely value near 
(3 = 2.2 (Schnee et al. 2010). Another example is the PRONAOS (PROgramme NAtional 
d' Observations Submillimertiques; Lamarre et al. 1994) balloon-borne experiment which 
finds that the opacity index has an inverse temperature dependence with (3 > 2 for T < 12.7 
K (Dupac et al. 2003). The lower bound of our opacity index overlaps with the PRONAOS 
opacity index curve (f3 = 1/(0.40 + 0.0079T); Dupac 2009) for the typical (T los ) found in our 
dust continuum models (see Figure 5). The PRONAOS results are not unique as an inverse 
temperature dependence of (3 and opacity indices greater than two at low temperatures has 
also been seen in far-infrared and submillimeter observations from the ARCHEOPS balloon- 
borne experiment (Desert et al. 2008). 
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The tendency of the opacity law toward (3 = 2 was noted from early submillimeter 
observations and is thought to originate from behavior of the complex dielectric function 
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e (Re(e)=const, Im(e)oc v\ Wickramasinghe 1967) of the grains at wavelengths far from 
resonances in the grain materials (e.g. Gezari et al. 1973). None of the popular opacity 
models used in protostellar dust continuum modeling or modeling of ISM dust predict power- 
law indices greater than two (i.e., OH2 (3 smm = 1.35, OH5 (3 smm = 1.85, OH8 (3 smm = 1.88, 
WD5.5 f3 smm = 1.69); see Shirley et al. 2005). While our individual submillimeter opacity 
constraints overlap the OH5 model at both 450 and 850 /xm due to the uncertainty in the 
physical model that best fits B335, the opacity index must be determined using the same 
physical model at both wavelengths and the resulting f3 smm is too steep to be consistent with 
the OH5 model. 

There is a class of amorphous silicate dust models which include phonon difference 
processes (Disordered Charge Distributions and localized Two Level Systems; see Schlomann 
1964 and Phillips 1987) that result in /3 submm > 2 at low temperatures (e.g., Meny et al. 
2007). These processes have been used to explain the anti-correlation between (3 smm and 
temperature observed by the PRONAOS experiment (Boudet et al. 2005). Our (3 smm range 
is consistent with the predicted opacity index from Boudet et al. (2005) and Meny et al. 
(2007) for the typical TJ os < 10 K derived from the dust models. Ultimately, our results 
should be tested by reproducing this analysis for observations with the new generation of 
submillimeter cameras (e.g., SPIRE, LABOCA, SCUBA2) and using better constraints from 
more sophisticated (e.g., multi-dimensional) dust continuum models of B335. 



5. Conclusions 

We have determined the opacity ratio from the slope of a plot of submillimeter intensity 
versus near-infrared color excess toward B335. The submillimeter intensity along each line-of- 
sight in the correlation is corrected for the non-isothermal temperature profile by the quantity 
1/P n which is related to the density- weighted average Planck function. We find opacity ratios 
of Msa = (3.21 - 4.80)1°;^ x 10" 4 and = (12.8 - 24.8)j£| x 10~ 4 for a ratio of total to 
selective K-band extinction of R K = 1.59 ± 0.12. The range in values corresponds to the 
uncertainty in the physical model for the envelope of B335. The submillimeter opacity power- 
law index is (3 smm = (2.18 - 2.58) IqJo- For an average 2.2 /iin opacity of 3800 ± 700 cm 2 g _1 , 
we find an opacity at 850 and 450 /im of k^o = (1.18 — 1.77)^Q24 and K450 = (4.72 — 9.13)tgg 8 
cm 2 g~ 1 . These opacities statistically agree with the popular theoretical ratios of Ossenkopf 
and Henning for coagulated ice grains with thin mantles (0.65 — 0.97)k^ 5 at 850 /xm; 
however, our derived opacity index ((3 smm ) is steeper than predicted by the OH5 model 
{Pews = 1-85). This comparison of near-infrared color excess and submillimeter emission 
probes the opacity on scales of 15" to 75", and does not find evidence for a large scale 
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variation in the opacity on those scales. We confirm a disagreement between the best fitted 
dust radiative transfer model and the best-fitted molecular line radiative transfer model. 
Improvements in the estimate of the opacity ratios and submillimeter opacities may be 
made with more sophisticated, milti-dimensional modeling of the dust continuum emission 
such as variable dust opacities in the inner envelope where desorption of CO and other 
molecules may change the optical constants of grains. The techniques used in this analysis 
should be applicable to far-infrared and submillimeter observations of B335 with the Herschel 
Space Observatory. With the commissioning of new, sensitive bolometer cameras, such as 
LABOCA and SCUBA-2, combined with observations with large format infrared CCDs on 
large aperture telescopes (JWST, Keck, etc.), it will be possible to extend this method to 
study the dust opacity ratio around other Class protostars. 
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Fig. 1. — Greyscale 850 /iin image of B335 with positions of background stars observed with 
NICMOS indicated by the small black crosses. All positions are in the epoch J2000.0. The 
850 (Ltm (black) contours starts at 2a (30 mJy/beam) and then are spaced at 10% of the 
peak (101 mJy/beam). The red and blue contours trace the outflow wings derived from CO 
J = 2 — > 1 observations and start at 1.5 K km/s and increase by 1.5 K km/s (Stutz et al. 
2008). The extent of the outflow cavity used to exclude background stars is shown by the 
two solid black lines. The SCUBA beam (lower left) and SMT CO beam (lower right) are 
displayed. 
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Fig. 2.— Updated best-fit dust model for B335, 2.4x Harvey BPL OH8. Panel (a) displays 
the scaled-Harvey broken power-law density profile (blue) and the resulting dust temperature 
profile (red). The fit to the SED is shown in panel (b). The histogram showing the location 
of background stars that are used to constrain the opacity ratio is shown in panel (c). 
The distribution is strongly peaked just below 10 4 AU. The fit to the submillimeter intensity 
profiles at 850 and 450 /zm are shown in the bottom panels. The red curves are the dust model 
profiles while the black curves are the beam profiles determined from Uranus observations 
bracketing the B335 observations. The intensity errorbars account for statistical uncertainty 
in the intensity as well as azimuthal variations in intensity within each annulus. 
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Fig. 3. — The 850 fim intensity plotted versus (H - K) color excess. P n is related the 
density-weighted Planck function. The solid line is the linear regression with the mean slope 
and intercept from the Posterior distributions. Histograms of the intercept (a) and slope (b) 
distributions are shown in the insets. «85o/ K 2.2 = b 850 /RK where we have assumed Rk = 1-59 
(see §2). 
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Fig. 4. — The 450 fxm intensity plotted versus (H - K) color excess. P n is related to the 
density-weighted Planck function. The solid line is the linear regression with the mean slope 
and intercept from the Posterior distributions. Histograms of the intercept (a) and slope (b) 
are shown in the insets. ^450/^2.2 = b^o/Rx where we have assumed Rk = 1-59 (see §2). 
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Fig. 5. — k 850 /k 2 .2 versus the average line-of-sight dust temperature. k 850 /k 2 .2 is determined 
from the slope of linear regressions (&850 see Figures 3 and 4) divided by Rk = 1-59 (see 
§2) for different physical models (n(r), T(r)). The solid line is for dust models with an 
isothermal envelope (T(r) = constant). The plotted points are for various dust models with 
calculated T(r). Symbols refer to: scaled Harvey broken-power law (triangles), scaled Shu 
model (squares), and Evans-Shu model (circle). Color indicates the dust opacity model used 
to calculate T(r): OH8 mo d (green), OH5 moc ; (blue), and OH2 mo d (red) (§3). All models 
except for the Evans-Shu model are scaled in density to match the observed 850 jttm flux. 
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Fig. 6. — Theoretical kappa ratios with the observed opacity ratios from this work (B335, 
shown in red), Bianchi et al. (B68, 2003), and Kramer et al.(IC 5146, 1.2mm 1998, 850 /im 
2003). The Kramer et al. 850 /im points have been shifted slightly in wavelength for clarity 
WD = Weingartner & Draine (2001) for Ry = 5.5. Mathis refers to the parametrization 
by Mathis (1990) of the ISM empirical dust model. OH model profiles are not included in 
this Figure since self-consistent scattering opacities were not determined by Ossenkopf & 
Henning (1994). 
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Fig. 7. — Our constraints on the submillimeter dust opacity assuming and absolute value 
of the opacity at 2.2/rni of 3800 ± 700 cm 2 g" 1 . WD = Weingartner & Draine (2001) for 
Ry = 5.5. OH = Ossenkopf & Henning (1994) for (2) no ice mantles, (5) thin ice mantles, 
and (8) thick ince mantles (N. B., the Young & Evans 2005 modification to the OH opacities 
does not affect the submillimeter opacities). Mathis referes to the parametrization by Mathis 
(1990) of the ISM empirical dust model. 
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0.26 


2.01 


43.20 


1.8 


8.52 


± 


0.73 


7.34 


± 


0.39 


1.18 ± 0.68 


36.5 


± 


2.3 


6.4 ± 4.3 


2.8 x Evans SHU OH8 d 


4.06 


12.37 


32.29 


2.7 


9.24 


± 


0.82 


6.00 


± 


0.32 


1.11 ± 0.56 


26.2 


± 


1.7 


5.8 ± 3.2 


11. Ox Harvey BPL WD5.5 


0.76 


1.54 


17.91 


2.7 


9.18 


± 


0.71 


7.61 


± 


0.37 


-0.75 ± 0.65 


39.2 


± 


2.1 


-1.0 ± 0.4 


3. lx Evans SHU OH5 


1.98 


10.42 


44.57 


2.1 


8.85 


± 


0.76 


6.73 


± 


0.36 


1.06 ± 0.64 


32.0 


± 


2.0 


5.1 ± 3.8 


1. 6 x Evans SHU OH2 


0.75 


6.00 


43.41 


1.5 


8.91 


± 


0.84 


6.18 


± 


0.33 


1.67 ± 0.60 


27.1 


± 


1.8 


10.8 ± 3.5 


14.0 x Evans SHU WD5.5 


0.23 


17.91 


37.45 


2.3 


9.35 


± 


0.72 


6.85 


± 


0.36 


-0.16 ± 0.62 


33.4 


± 


1.9 


-5.2 ± 3.3 


1.0 x Evans SHU OH5 


7.72 


15.12 


53.49 


1.9 


9.74 


± 


0.93 


5.11 


± 


0.27 


1.33 ± 0.47 


20.4 


± 


1.4 


7.@g± 2.7 



t 



a Harvey BPL = Harvey et al. (2003b) broken power-law n(r) and Evans SHU = Evans et al. (2005) Shu-infall n(r). The numbers 
refer to scaling factors multiplied into the density or the opacity. All models except 1.0 x Evans SHU OH5 are scaled to match the 
observed flux at 850 /zm. 

b Luminosity integrated from A > 60 /im. The B335 SED is published in Shirley et al. (2002) plus additional points from Spitzer 
Space Telescope observations at 70 /im (S = 15.4 ±2.1 Jy in a 70" aperture) and 160 /im (S = 68.7 ± 15.6 Jy in a 100" aperture) 
from Stutz et al. (2008). 

c Best-fitted dust continuum model. 

d Best-fitted molecular model that matches flux at 850 //m. 

e a and b refer the the intercept and slope of the linear regression at 850 and 450 fxm. 



